#include "parafluidsmoothingkernel.h"
#include "parafluidsph.h"

#include <cmath>

using namespace ParaFluid;

const float Pi = 3.14159265;

float SmoothingKernel::valuePoly6(const float & h, const float & r)
{
	if(0 <= r && r < h)
	{
		float factor = 315.0 / (64.0 * Pi * h * h * h * h * h * h * h * h * h);
		float temp = (h * h) - (r * r);
		return factor * temp * temp * temp;
	}
	else
	{
		return 0;
	}
}

Vector SmoothingKernel::gradientSpiky(const float & h, const Vector & R)
{
	float r = R.length();
	if(0 < r && r < h)
	{
		float factor = 45.0 / (Pi * h * h * h * h * h * h);
		return R * factor * ((((h * h) + (r * r)) / r) - (2.0 * h));
	}
	else
	{
		return Vector(0.0, 0.0, 0.0);
	}
}

float SmoothingKernel::laplacianViscosity(const float & h, const Vector & R)
{
	float r = R.length();
	if(0 < r && r < h)
	{
		float factorOfMagicalReduction = pow((h - r) / h, 3.0);
		return (45.0 / (Pi * h * h * h * h * h * h)) * (h - r) * factorOfMagicalReduction;
	}
	else
	{
		return 0.0;
	}
}
